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INTRODUCTION 

The  verification  of  numerical  schemes  for  solving  the  equations  of  inviscid  and  viscous  compressible  unsteady  flow  equations  is 
limited  to  a small  number  of  analytical  solutions  of  the  equations  governing  the  one-dimensional  unsteady  flow  including 
moving  discontinuities.  Among  them  the  most  important  were  given  first  by  B.  Riemann  (1859-1860)  and  later  by  W.J.M. 
Rankine  (1870),  P.H.  Hugoniot  (1887),  Lord  Rayleigh  (1910)  and  G.I.  Taylor  (1910).  The  scope  of  the  present  chapter  is  to 
overview  the  analytical  solutions,  serving  as  test  case  for  the  accuracy  of  the  numerical  schemes.  It  is  worth  noting  that  the 
analytical  solutions  are  of  importance  for  Euler  and  Navier-Stokes  equations  for  laminar  flow  and  does  not  give  any  indication 
for  the  behaviour  of  the  numerical  schemes  in  the  prediction  of  turbulent  flows.  For  each  analytical  solution  a corresponding 
FORTRAN  program  is  attached. 

LINEAR  ADVECTION  EQUATION 

The  linear  advection  equation  is  used  as  a simple  model  for  contact  discontinuities  in  fluid  dynamics,  and  is  the  simplest  model 
equation  for  the  representation  of  wave  propagation.  The  linear  advection  equation  is: 

u,  + c0ux  = 0 (1) 


where  c0  is  a positive  constant  (c0  > 0 ) called  the  velocity  of  the  wave.  The  general  solution  of  equation  (1)  is 

u(x,t)  = f(x-c0t)  (2) 

f(x)  is  an  arbitrary  function  defined  by  the  initial  conditions  of  the  problem 

u(x,0)=  f(x)  (3) 

The  solution  (3)  apparently  describes  a wave  motion  to  the  positive  x-axis,  since  the  initial  profile  f(x)  is  translated  unchanged 
in  shape  a distance  c0t  to  the  right  at  time  t [8]. 

Oscillatory  solution  of  linear  advection  equation  with  a discontinuity  in  the  derivative 
(case  LADV-1) 


The  problem  of  oscillatory  solution  of  linear  advection  equation  is  approached  using  the  initial  conditions  [8]: 


u(x,0)  = f(x;tx,k)  = 


-asin(knx),x  < 0 
x , x>0 


with  a = 0. 1 and  k = 6 


(4) 


At  x=-l  the  boundary  condition 


u(-l,t)  = 0 (4a) 

is  imposed.  The  solution  in  accordance  with  the  above  consideration,  is  presented  in  Figure  1 at  time  instances  t=0  and  t=0.5 

This  test  case  gives  information  about  the  capability  of  a computational  method  to  capture  an  oscillating  solution  with  a 
discontinuity  in  the  derivative. 
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Figure  1:  Oscillatory  solution  of  linear  advection  equation  with  a discontinuity  in  the  derivative 


Simulation  of  discontinuities  (case  LADV-2) 


To  simulate  discontinuities  we  again  solve  the  linear  advection 

u(x,0)  = f(x;uL,uR ) = |Ul’x <0  with 
uR,x>0 


equation,  this  time  using  piecewice  continuous  initial  data: 
uL  = 1 and  u R — 0 


(5) 


At  x=-l  the  boundary  condition 


u(-l,t)  = 0 (5a) 

is  imposed.  The  resultant  solution  is  a square  wave  travelling  with  speed  1 to  the  rigid,  as  it  is  depicted  in  Figure  2 . 
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Figure  2:  Contact  discontinuity,  advection  equation. 


BURGERS’  EQUATION 

Inviscid  Burgers’  equation. 

I'he  non  linear  first  order  equation: 
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u,  + uux  = 0 


with  Initial  Conditions 


u(x,0)  = f(x)  for  -oo  < x < oo 


(6) 


(6a) 


is  the  x-momentum  equation  without  pressure  gradient  or  other  external  forces  and  it  is  so  called  inviscid  Burgers’  equation. 

The  general  solution  of  the  above  equation  is  given  by: 

dx 

du  = 0 along  the  characteristic  - — = u . 

dt  v / 

expressing  that  u remains  constant  on  the  characteristic. 

For  the  above  initial  distribution  of  u(x,0)=f(x)  the  general  solution  is  concluded  in  an  implicit  form,  as  follows  [4]: 

u(x,t)  = f(x-ut)  (8) 

The  characteristics  have  a slope  proportional  to  l/f(x0)  in  the  (x,t)  plane,  where  Xo  is  a position  at  initial  state,  and  if  f'(x0)  is 
positive  , which  is  typical  for  an  expansion  profile,  they  will  never  intersect.  On  the  other  hand  for  a decreasing  initial 
distribution  of  u,  that  means  f '(x0)  < 0 , the  characteristics  will  intersect  as  for  a typical  compression  profile.  An  initial  profile 
with  decreasing  intensities  will  lead  to  a breakdown  of  a continuous  solution  and  to  the  appearance  of  a shock  discontinuity.  The 
shock  will  appear  at  the  time  instance  ts , when  the  tangent  to  u(x)  profiles  becomes  vertical: 


’ max|f'(x0)|  (9) 

The  shock  wave  velocity , us , satisfies  the  Rankine-Hugoniot  relations  and  is  equal  to: 


«s=^Uj  + Ul)  (10) 

where  ui,  u2  are  the  values  upstream  and  downstream  of  the  shock. 

At  this  point  the  following  three  types  of  initial  conditions  are  proposed  for  the  analyses  of  non  oscillatory  shock  capturing 
methods. 


Initial  shock  discontinuity  (case  IB-1) 

For  this  case,  the  Riemann  problem  for  Burgers’  equation  is  solved.  The  test  case  is  provided  by  an  initial  discontinuous 
distribution  : 


u(x,0)  = f(x;uL,uR)  = 


{uL,x  < 0 
uK,x>0 


with  uL  > uR 


(ID 


At  the  left  boundary  the  condition 


u(xL,t)  = uL  (11a) 

is  imposed.  The  solution  of  Burgers’  equation  gives  a shock  propagating  at  speed  (uL  + uR)/2  with  unmodified  intensity 
[u]=  uL  -uR  , as  shown  in  Figure  3.  If  uK  = -uL , the  shock  is  stationary  and  it  is  used  as  a non-linear  test  case  for  steady-state 
methods. 
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Figure  3 : Burgers’  solution  for  a propagating  discontinuity. 

Initial  linear  discontinuity  (case  IB-2) 

A different  initial  distribution  with  f ’(x0)  < 0 leads  to  the  same  shock  structure.  The  initial  linear  distribution  is: 


u(x,0)  = f(x;uL,uR)  = 


uL,  x<0 

X i X 

1-— l + uR— , 0<x<L  with  uL  > uR 
u- , x > L 


and  at  the  left  boundary  the  condition 


u(xL,t)  = uL 

is  imposed.  A shock  is  formed  at  time  instance 


‘s  = 


L 


Ul-Ur 


and  at  position  xs  = tsuL  = L + tsuR  . The  solution  of  Figure  4 for  t > ts  is  : 


u(x,t) 


uLforx 
uR  forx 


< 


> 


u,  +UR 


t 

t 


(12) 


(12a) 


(13) 


(14) 


Burgers’  equation  for  a rarefaction  wave  (case  IB-3) 

Burgers’  equation  with  the  following  initial  conditions  gives  a propagating  rarefaction  wave. 

f u,  x < 0 

u(x, t = 0)  = f (x;uL , uR ) = ( with  uL  <uR 

uH  x>0 
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(15) 


Between  points  uLt  < x < uRt  the  solution  is  not  determined  by  the  intersection  of  characteristics.  So,  a continuous  solution  is 
possible  in  the  following  form  (Figure  5): 


u(x,t) 


x/t  < uL 
uL  <x/t<uR 
x/t>u8 


(16) 


Figure  5 : Initial  state  and  time  evolution  of  a propagating  rarefaction  wave  as  a solution  of  Burgers'  equation. 


Viscous  Burgers’  equation  (case  VB-1) 


The  complete  nonlinear  Burgers’  equation  adding  a viscous  term  is  : 


u,  + uu,  = vu„ 


with  Initial  Conditions 


u(x,0)  = f(x)  for  -co<x<oo 


(17) 


(17a) 


The  “Viscous”  Burgers’  equation  serves  as  a model  equation  for  the  boundary-layer  equation,  the  “parabolized”  Navier-Stokes 
equations  and  the  complete  Navier-Stokes  equations. 

The  problem  with  the  following  initial  values: 


u(x,0)  = f(x;uL,uR)  = 


uL,x  < 0 
uR,x>0 


with  uL  > uR 


(18) 


and  boundary  conditions 


u(x  = -co,  t)  = uL , u(x  = oo,  t)  = u. 


(18a) 
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has  a solution  of  the  form: 


U 


u,  + u„ 


1 + hexp 


2v 


(x-Ut) 


where: 


je^dC 

-<x-uRt)//4vt 

]e-^dC 

(x-uLt)/y4vt 


(19) 


The  diffusing  shock  still  propagates  with  the  “inviscid”  velocity  equal  to  U.  Due  to  viscosity  effects  the  inviscid  discontinuities 
are  transformed  into  continuous  shaped  “steps”,  as  it  is  shown  for  the  test  case  of  Figure  6. 


Figure  6 : Shock  wave  solution  of  viscous  unsteady  Burgers’  equation  (UL  = 2 , UR  =1). 


UNSTEADY  EULER  EQUATIONS 

Reflection  of  a moving  shock  on  a closed  boundary  (case  RMS-1) 

The  general  discontinuity  equations  for  the  moving  shock  are  shown  in  the  literature  reference  [1]  to  be: 


P =_u 
P “ 


- = 1 


2 

7 + 1 


i \ 
T 

y 


i=i +-?l 

P 7 + 1 


vc 


\ 

) 


where  u,  u denotes  the  relative  velocities  in  front  and  behind  the  moving  shock 


(20) 

(21) 


(22) 


u = o - w, 


u = u - w 


(23) 
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while  V),  U are  the  absolute  velocities  and  w the  velocity  of  the  shock  front.  The  velocity  of  sound  c in  front  of  the  moving 
shock  wave  is  given  by: 


c = 


(24) 


The  test  problem  proposed  here  is  the  reflection  of  a shock  wave  moving  with  constant  velocity  towards  the  closed  boundary  of 
a tube  (figure  7).  The  fluid  behind  the  shock  wave  moving  to  the  left,  has  a velocity  with  the  absolute  value  u0  (o  = -o0 , 


o0  > 0 ),  pressure  and  density  p„,p0  ( p = p0,P  = Po ) respectively,  so  that  the  velocity  of  sound  be 


. The  given  data 


are  p0,p0,u0,  which  are  the  initial  conditions  of  the  problem  at  starting  point.  The  reflection  condition,  which  has  to  be 


satisfied  is  u = 0 . 


Figure  7 : Reflection  of  a moving  on  the  left  with  constant  velocity  shock  wave  on  the  left  closed  boundary  of  a tube. 

The  velocity  of  the  shock  front  after  the  reflection,  the  pressure  and  the  density  behind  the  reflected  shock  wave,  are  given  by  the 
following  relationships  [2]: 


(25) 


P = Po  1 


27 

Y + l 


(u„+w)2 


YPo 


Po"1 


(26) 


- o„  + w 
p = — Po 


(27) 


Analytical  solutions  for  the  unsteady  inviscid,  non-conducting  fluid  conservation  equations 

By  neglecting  viscous,  heat-conduction  effects  and  field  forces,  the  unsteady  compressible  conservation  equations  of  mass, 
momentum  and  energy,  in  one  dimensional  conservative  form,  which  will  be  referred  to  as  Euler  equations,  have  the  following 
form: 


Continuity  equation: 


dp  , d(pu)  _ Q 
dt  dx 


Momentum  equation: 


(28) 


d(pu)  | d(pu2  +p)_Q 
dt  dx 


(29) 
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Energy  equation: 


1 - 
p|  e + — u' 
1 2 


1 2 P 
pu  e + — u + — 

' 2 p 


= 0 


<3t  3x 

The  above  system  is  closed  by  the  constitutive  equation,  in  the  form  of  Gibbs  relation , namely: 

Tds  = de  + pd| 


(30) 


Which  leads  to  the  entropy  relation: 


5(ps)  d(pus) 


0 


dt  dx 

In  the  space  - time  plane  the  transformation  of  independent  (x,t)to  the  new  variables  (i;,  r|)  (Figure  8)  is  introduced: 


(31) 


(32) 


Figure  8:  Sketch  of  the  two  families  of  characteristics  in  the  (x,t)  plane. 

£>  = ^(x,t),Ti  = t|(x,t)  (33) 


so  that: 

, dx  dx 

c = const:  — =u-c,  n = const:  — = u + c 
dt  dt 


where  c is  the  isentropic  velocity  of  sound: 


Then  it  can  be  shown  that  the  system  of  governing  equations  in  terms  of  (^,r|)  takes  the  form: 


dx 

n = const : — = u + c,  u + to  = const 
dt 

. dx 

t = const : — = u-c,  u-co  = const 
dt 


where: 


0)  = 


(34) 


(35) 


(36) 


(37) 
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The  E,  = const.,  q = const,  are  the  two  families  of  characteristics  which  are  wave  fronts  of  the  kinematic  discontinuities,  lhe 
kinematic  discontinuities  in  the  considered  one-dimensional  case  correspond  to  lines  across  which  the  first  derivative  of  the  flow 
quantities  are  discontinuous  , while  the  flow  quantities  are  continuous. 

The  last  equation  for  ideal  gas  of  constant  coefficients  of  specific  heat  of  ratio  y reduces  to  the  form: 


co  = 


2 

c 

7-1 


(38) 


The  reduced  system  of  equations  can  be  solved  analytically  in  certain  problems  such  as  moving  piston  case  (two  sub-cases: 
expansion  and  compression)  and  the  so  called  Riemann  or  shock  tube  problem  which  includes  a shock  wave,  a contact 
discontinuity  and  an  expansion  wave  [1],  [3]. 

(a)  Expansion  flow  behind  a moving  piston  (case  MP-1) 

A piston  is  considered  to  move  towards  the  negative  x-direction,  figure  9.  This  results  to  an  expansion  of  the  gas  behind  the 
piston. 


Figure  9:  Expansion  flow  on  the  right  of  a left  moving  piston 

The  flow  is  studied  in  the  plane  x-t  where  the  path  of  the  piston  lA=-TJis  shown.  As  an  initial  value  of  the  problem  a gas  with 
zero  velocity  and  a constant  temperature  (constant  velocity  of  sound)  will  be  considered. 

u(x,0)  = 0,  c(x,0)  = c0  x > 0 (39) 

As  a boundary  condition  the  equality  of  the  gas  velocity  to  the  velocity  of  the  piston  is  taken  into  account: 

u(xP-t)=up(t)=ltr  (40) 

Starting  from  the  initial  values  c = c0,u  = 0 (for  t = 0 ) we  could  easily  observe  that  in  the  open  region  between  the  positive  x- 
axis  and  the  positively  inclined  characteristic  ( q = const.  ) originating  from  point  1 both  families  of  characteristics  are  straight 
lines  since  they  originate  from  the  positive  part  of  the  x-axis  where  the  constant  initial  conditions  are  valid. 

The  movement  of  the  piston  affects  the  flow  field  left  of  the  characteristic  originating  from  point  1.  Since  all  the 
q = const  characteristics  originate  from  the  positive  x-axis  where  the  initial  conditions  are  valid,  across  all  the  negatively 
inclined  characteristics  the  following  relation  is  valid: 

2 2 

£ = const,  u c = cn 

Y - 1 y-1  (41) 


All  the  E,  = const  characteristics  end  on  the  piston  path,  so  that 
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£ = const  u„ — — c = — c 

p y-1  p y-1 

while  for  the  q = const  characteristics  the  following  relation  is  valid: 

2 2 


(42) 


n = const  uh c=u„+ c„ 

y-1  p y-1  p 


(43) 


By  substracting  the  previous  two  equations  for  the  region  between  the  piston  path  and  the  positively  inclined  characteristic 
originating  from  point  1 , the  following  relation  along  the  positively  inclined  characteristics  is  valid  : 


2 2 

q = const  u + c = 2u  + c0 

y-1  y-1 


(44) 


From  the  relations  (42)  and  (43)  wre  conclude  that: 

q = const,  u = up,c  = cp  (45) 

Thus  the  values  of  u and  c on  each  r|  = const  keep  constant  and  this  family  of  characteristics  are  straight  lines  with  the 
following  inclination  : 

dx  y + 1 

11  = const-  — = u,+cp=c0+— up  (46) 


Note:  From  equation  results  a limit  maximum  piston  velocity. 


K)m„  =- 


y-1 


cp=0 


(47) 


With  a piston  velocity  increasing  further,  a cavitation  zone  is  developed  behind  the  moving  piston  (where  pressure  vanishes), 
The  cavitation  area  is  located  in  the  region  between  the  piston  path  and  the  positively  inclined  characteristic  with  inclination 
w : 


dx  y + 1 

w = — = c0  + 

dt  0 2 


(u  ) 

' P 'mi 


(48) 


This  limit  velocity  w » -5c0  is  quite  higher  than  the  maximum  isentropic  steady  flow  velocity  in  vacuum  ( umK  = -J~5c0 ).  This 

shows  the  basic  differences  between  the  steady  and  unsteady  flows.  Of  course,  we  should  note  that  these  regions  are  on  the  limit 
of  continuums  mechanics  validity. 

The  whole  phenomenon  can  be  considered  as  simple  wave.  This  is  the  case  of  the  wave  when  one  family  of  the  characteristics 
are  straight  lines. 

Special  case  I: 

In  the  special  case  of  a piston  moving  with  constant  velocity,  so  that  the  piston  path  is  a straight  line.  The  expansion  region  is 
limited  to  an  angle  around  the  axis  origin,  so  that  all  the  gas  particles  hold  the  constant  velocity  of  the  piston.  The  wave  is  called 
central  simple  expansion  wave. 
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Figure  10  : Central  expansion  wave  behind  a piston  moving  with  constant  velocity 


Special  case  II: 

For  the  special  case  of  the  constant  accelerating  piston,  for  which  the  velocity  argument  is  linearly  increasing  with  time  (U  and 
t0  are  constants): 


(49) 


leads  to  a folly  analytical  expression  of  the  flow  velocity  and  the  isentropic  velocity  of  sound  in  the  expansion  area,  which  are 
given  by  the  relations: 


u(x,t)  = — 

y 


2 t„ 


yl  (c° 


2 t„ 


2 1 3 

2 U / , 

+- — (x-e0t) 

y 


(50) 


c(x,t)  = 


Y -1 

co+~ u 


(51) 


The  distribution  of  velocity  as  a function  of  the  space  variable  x is  shown  in  figure  1 1 for  various  time  levels. 
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Figure  11:  Velocity  distribution  as  function  of  x for  various  times  for  the  case  of  gas  expansion  due 
to  the  movement  of  a piston  with  a constant  acceleration  ( U = 100m/s,t0  = ls,c0  = 400m/s,y  = 1.4  ). 


(b)  Compression  flow  in  front  of  a piston  moving  in  a non-moving  gas  (MP-2) 

The  theory  of  isentropic  flow  for  the  compression  flow  in  front  of  a piston  obeys  the  same  analysis  as  the  expansion  one  (eqs 
(38)-(47)).  One  should  remark  in  this  case  that  the  t|  = const  characteristics  converge  and  form  an  envelope  (Figure  12).  The 
envelope  can  be  in  general  shown  that  appears  at  earliest  time  at  the  point  (xc  ,tc) , that  is  defined  as: 


2Cq 

(y  + l)u  (0)’ 


t„  = 


2c„ 


(y  + l)u  (0) 


(52) 


2 c _ 2 


co 


Figure  12  : Compression  flow  in  front  of  a piston  moving  in  a non-moving  gas 
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After  this  point  and  for  ( t > tc ) a moving  shock  wave  appears,  which  propagates  in  the  gas  at  rest  in  the  same  direction  with  the 
piston  and  the  flow  is  anisentropic.  An  analytic  description  of  the  flow  field  succeeds  for  the  following  two  cases: 

Special  case  III 

In  the  case  that  the  piston  moves  with  constant  velocity  in  a gas  at  rest: 

up  = U , U = const  (53) 

so  that  the  piston  path  is  a straight  line,  figure  13,  the  p = const  characteristics  are  parallel  to  each  other  and  a shock  forms, 
which  propagates  with  constant  velocity  us  higher  than  that  of  the  piston: 


dxs  _ y + 1 

~dT~~ 


U + 


+ c 


2 

0 


(54) 


Special  case  IV 


For  the  case  of  the  compression  flow  of  a gas  in  front  of  a moving  piston  with  a velocity  linearly  increasing  with  time  (constant 
accelerating  piston): 


(55) 


(U,  t0  constants),  we  lead  as  in  the  case  II  to  a fully  analytic  expression  for  the  flow  field  velocity  in  the  isentropic  region,  that  is 
before  the  appearance  of  the  shock  wave.  The  expression  of  the  velocity  in  the  isentropic  compression  region,  can  be  shown  to 
be  similar  to  the  expression  of  the  expansion  flow: 


2 U, 
— — (x- 
7 to 


c(x,t)=  c0 


and  the  point  (xc,tc)  is  calculated  as: 


x _ 2coto  _ 2c0t0 

c (y  + l)u  ’ ' (y  + l)u 


(56) 


(57) 


(58) 


The  distribution  of  velocity  as  function  of  the  space  variable  x for  various  times  t is  shown  in  figure  14. 


22 


Figure  14  : Distribution  of  the  velocity  as  a function  of  x for  various  t in  the  isentropic  compression  area 
of  a gas  in  front  of  a constant  accelerating  piston  ( U = 1 00  m / s,  t0  = 1 s,c0  = 400  m / s,y  =1.4). 


(c)  Moving  shock  in  a shock  tube  or  Riemann-problem  (case  ST-1) 

The  shock  tube  problem,  often  called  Riemann  problem,  is  the  flow  owing  to  the  abrupt  removal  of  the  valve  which 
separates  a high  pressure  gas  from  a low  pressure  gas  in  the  shock  tube  [1],  [6].  The  resulting  wave  effect  of  the  propagation  of 
the  discontinuity  and  the  relating  nomenclature  are  shown  in  figure  15.  The  diaphragm  located  in  position  x = 0 at  time  t = 0 
separates  the  space  of  high  pressure  pH  s p4  from  the  space  of  lowpressure  pL  s p,  . Thus,  the  basic  parameter  of  the  flow 
effect  is  the  pressure  ratio: 


Ph  _ P< 
Pl  Pi 


(59) 


The  two  parts  of  the  tube  may  have  also  different  temperatures  ( T4 , T,  ) and  different  gases  ( R4 , R , ). 

At  initial  time  the  pressure  distribution  is  a step  distribution,  figure  15a.  This  causes  the  separation  of  the  problem  in 
two  problems,  the  propagation  of  a shock  wave  in  the  low  pressure  gas  pL  and  the  propagation  of  an  expansion  wave  in  the 


high  pressure  gas  pH  . The  state  behind  the  moving  shock  wave  is  indicated  by  the  index  2 and  the  state  behind  the  expansion 
wave  with  the  index  3.  The  interface  between  the  states  2 and  3 is  a contact  discontinuity.  This  is  the  contact  point  of  the  two 
gases,  initially  separated  by  the  diaphragm,  and  have  different  temperatures  and  densities.  On  the  other  hand  they  should  have 
the  same  velocities  and  pressures.  The  basic  problem  is  how  the  flow  quantities  can  be  calculated  with  a given  initial  pressure 
ratio. 


Introducing  the  following  expressions: 


Pl 


.1+1 

Y-l 


(60) 


we  have  firstly  the  relations  connecting  the  quantities  on  both  sides  of  the  moving  shock  wave: 
Moving  shock  wave  relations: 


p2  1 + aP  , TI  . . , 

= ( Hugomot  relation) 

Pl  a + P 

(61) 

u2-o,.  f 2 V P-1 

(62) 

^(r-Oj  (1  + ap)i 
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c,  fp“  + Pf 

cL  V l + aPy 


W-CL 


h aP 


(63) 


(64) 


The  pressure  and  velocity  on  both  sides  of  the  contact  discontinuity  are  equal. 


diaphragm 


y///////////////////////////////////////////w^ 


H High  pressure 


Low  pressure 


(1)  L 


(4) 

V777777777777777777777Z?77777Z7777?7w7777f777W!7777777777777777777777777777777777Z777777f7Z 
(a)  t=0 


Equations  of  the  contact  discontinuity: 


Pr  P2 

U3  =Oj  = V 


(65) 


v is  the  velocity  of  the  contact  discontinuity. 
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Pressure  and  velocity  in  the  regions  3 and  H = 4 are  related  with  Riemann  conditions  on  positive  inclined  characteristics 
(r|  = const) . 

Riemann  equations  across  positive  inclined  characteristics: 


After  eliminating  unknowns  from  equations  (61) , (64),  (65)  the  following  equation,  which  has  as  unknown  the  pressure  relation 
P,  is  obtained: 


P-1 
(l  + aP)i 


T-l  " 

CH  2 

l_f Pkp 

O 

t- 

1 

tT 

X 

The  solution  of  the  implicit  algebraic  equation  (66)  is  accomplished  by  numerical  integration. 


(67) 


Example: 

For  the  following  values  of  the  variables: 

Ph=105,ph=1,  t>H  — 0,  PL=10\  pu  =0.125,  oL=0,y  = 1.4 

the  solution  of  equation  (66)  gives  the  following  values  of  the  unknown  parameters: 

P = 3.0313,  p2  = 30313,  o2  = v = 203,  w=544 

In  figure  16  the  variation  of  flow  quantities  are  shown  for  the  time  instant  t = 6,1 10~3  . 


25 


0.00  2.00  4.00  6.00  1.00  10.00 

x - coordinate  (m) 


Figure  16  : Distribution  of  flow  quantities  at  t=6,l  ms  for  pH  = 105 , pH  = 1,  oH=0,  pL  = 104,pL  = 0.125,  uL=0,y  = 1.4 


ANALYTICAL  SOLUTIONS  FOR  THE  UNSTEADY  COMPRESSIBLE  LAMINAR 
FLOW  FOR  HEAT  CONDUCTING  FLUID 

The  existing  analytical  solutions  for  the  unsteady,  compressible  laminar  flow  for  a heat  conducting  fluid  concern  only  the 
Lighthill’s  approximation  of  finite  amplitude  sound  waves,  for  the  case  of  one-dimension.  For  the  derivation  of  this  theory  we 
refer  to  the  review  of  M.J.Lighthill  [5].  This  theory  leads  to  the  following  equation  for  the  velocity  u derived  from  the  equations 
of  continuity,  momentum  and  energy  after  simplifications  coming  from  the  assumption  of  sound  waves  of  finite  amplitude  for  a 
perfect  gas  of  constant  y : 


Su  ( 

aT° 


5u 

5x 


5£>-u 
2 3x2 


(68) 


5 is  the  Lighthill’s  diffusion  coefficient  for  the  sound  waves  of  finite  amplitudes  propagating  to  the  positive  - x direction: 


s_  do  t 

Po  PoCp 


(69) 


u„  = 2u  i u' , where  p p'  are  the  dynamic  and  volumetric  viscosities  of  the  gas  , X0  is  the  coefficient  of  heat  conductivity, 
cp  is  the  heat  coefficient  for  constant  pressure  and  p0  is  the  undisturbed  density. 


Through  the  transformation: 


u„=0,  c = c0 


_ V + 1 

X = x-c0t,  u = c0+J — -u 
4 


(70) 


(71) 
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the  differential  equation  can  be  transformed  to  the  following  non-linear  Burgers’  equation  (1940),  which  is  the  simplest  non 
linear  equation  describing  convective  effects  combined  with  diffusive  one  : 


3u  _ 9u  _ 8 8Ju 

aT  “ax-yax7 


(72) 


The  general  solution  of  the  Burgers’  equation  is  defined  by  the  introduction  of  a new  dependent  variable,  the  function  <|> 
defined  by  the  following  relations,  satisfying  the  Burgers’  equation: 


cty  _ a<(i  1 _2  8 5u 

— — U,  — 11 /"'I'}  \ 

ax  a 2 2 ax  (7J) 


The  differential  equation  for  <j)  leads  by  elimination  of  u from  the  previous  equations: 

aj>  = i f a<t>  V | 5 a2<t> 
at  ~ 2 1, ax J +2  ax2 


(74) 


By  introducing  again  a new  function  vp  : 


4>  =Slnip 


(75) 


the  equation  which  satisfies  vp  is  the  standard  linear  equation  for  heat  transfer  by  conduction: 

avp  _ 1 5 a2tp 

a ~ 2 ax2 


(76) 


So,  all  the  known  solutions  of  the  linear  heat  conduction  equation  are  at  the  same  time  solution  of  the  Burgers’  equation,  of 
course  in  the  transformed  variables. 

The  relations  connecting  vp  and  u are: 


u = — 


5 dvp 

vax’ 


' 1 00 

ip  = exp  — JudX 

x 


(77) 


As  a first  example  we  refer  to  the  initial  value  problem  owing  to  Laplace.  When  the  initial  value  of  the  wave  form  is  given  by 
u(x,0) , then  the  solution  is  defined  by  the  integral: 


tf(X,  t)  = -^==  Jvp(Y,0)exp 


(X  Yf 
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dY 


(78) 


u(x,t)= 


]- 

r'nf 

>(Y,0)dY-<X-fr] 

-dY 

j 

K| 

[ Ju(Y,0)dY  _ iX  - vf  | 

[i  2t  j 

(79) 


INDICATIONS  FOR  USE  OF  THE  EXAMINED  TEST  CASES 

The  above  test  cases  are  frequently  used  for  testing  various  properties  of  the  examined  numerical  schemes.  Specifically: 

• The  test  case  of  linear  advection  equation  (LADV-1)  that  results  in  oscillatory  solution  and  contains  discontinuity  in  the 
derivative,  can  be  used  to  test  the  diffusion  and  dispersion  properties  of  schemes  and  to  define  the  accuracy  of  the  scheme  on 
smooth  functions  of  the  wave  number  k. 

• The  linear  advection  equation  of  propagating  discontinuity  with  the  velocity  a (LADV-2)  is  important  with  regard  to 
properties  of  the  schemes  at  handling  propagating  discontinuities.  If  the  discontinuity  is  an  expansion  shock,  the  numerical 
scheme  can  propagate  and  dump  the  expansion  through  an  introduced  entropy  condition  or  any  other  form  of  dissipative 
mechanism. 
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• The  behaviour  of  invisvid  Burgers’  equation  against  non  linearities  is  representative  of  the  exantination  of  Euler  equations 
behaviour  due  to  the  non  linear  term  of  Burgers’  equation.  Inviscid  Burgers’  equation  (IB-1)  with  initial  shock  discontinuity 
gives  information  about  the  capability  of  the  numerical  scheme  on  shock  capturing  with  the  correct  shock  propagating  speed 
(time  accurate  scheme).  The  shock  capturing  without  the  presentation  of  non  physical  oscillations  in  the  vicinity  of  the 
discontinuity  ensures  the  monotonicity  of  the  numerical  scheme  or  its  property  for  total  variations  diminishing.  . Inviscid 
Burgers’  equation  (IB-2)  with  initial  linear  discontinuity  gives  information  about  the  diffusion  and  dispersion  properties  of  the 
numerical  scheme.  EB-2,  also  ensures,  as  IB-1  for  the  shock  capturing  capability  of  the  scheme  and  its  characteristics  about 
monotonicity.  Inviscid  Burgers’  equation  (EB-3)  for  a rarefaction  wave  is  used  to  test  the  additional  entropy  condition  that  is 
imposed  on  the  numerical  schemes  in  order  to  capture  expansion  shocks  for  inviscid  flow  equation.  Viscous  Burgers’  equation 
(VB-3)  serves  as  a model  equation  for  the  boundary-layer  equation  or  the  “parabolized”  Navier-Stokes  equations  and  examines 
the  capturing  of  diffusing  shock. 

• The  problem  of  shock  tube  presents  an  exact  solution  to  the  full  system  of  one-dimensional  Euler  equations  containing 
simultaneously  a shock  wave,  a contact  discontinuity  and  an  expansion  fan.  Consequently,  it  can  be  used  for  the  testing  of  all 
the  above  properties  of  numerical  schemes. 


SOFTWARE 

The  analytical  solutions  for  the  above  cases  are  also  presented  in  the  corresponding  FORTRAN  90  programs  that  have  been 
attached  to  the  present  paper.  These  programs  are: 


FORTRAN  program 

Input  Data  file 

Test  case 

advection.for 

advect.  int 

Linear  advection  equation 

inv Burgers. for 

inv  burg,  int 

Inviscid  Burgers’  equation 

vis  Burgers. for 

visburg.int 

Viscous  Burgers’  equation 

expansion.for 

expansion,  int 

Expansion  flow  behind  a moving  piston 

compresionfar 

compression,  int 

Compression  flow  in  front  of  a moving 
piston  in  a non-moving  gas 

shock  tube  for 

shock  int 

Moving  shock  in  a shock  tube  or  Riemann 
problem 
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